capture mkdir "${base_dir}/results/rd_plots"


*local bw 180
local folder rd_plots
local fe_opts rte // rte // rte // rtebyves // rte_ves 
local route_opts all_exits all_coast // all_exits //  // all_exits // all_coast      // socal_exits socal_coast   socal_exits  other_coast other_exits 
local vestype "Container"
local controls_resid //  IFO LS

local opts1 msymbol(O) 
local opts2 msymbol(D)

* recovering eca dates
levelsof eca_ind if eca_ind>0 , local(eca_inds)
foreach i of local eca_inds {
	qui sum date1 if eca_ind==`i'
	scalar eca`i'd = r(min)
	scalar eca`i'qy = yq( year(r(min)),quarter(r(min)) ) - yq(2009,1) 
	scalar eca`i'my = ym( year(r(min)),month(r(min)) ) - ym(2009,1) 
}

scalar list _all


local bin_var my
local omit_bin 592
local xline_opts = "`=eca1my' `=eca2my+1' , lpattern(dash) lcolor(black)" // `=eca3my' `=eca5my-1'  
* adjusting because of omitted bin

display "`xline_opts'"

/*
local bin_var qy
local omit_bin 197
local xline_opts = "`=eca1qy' `=eca2qy+1' , lpattern(dash) lcolor(black)" // `=eca3qy' `=eca5qy-1' 
*/


foreach rte_opt of local route_opts {
	
	* default options
	*local yvars d_km_below30 // SBind d_km_below24 d_km_above40 // ldist_eca2011 ldist_eca2009 ldist lf_cons // ldist_eca2011 ldist_eca2009 ldist lf_cons lcost_comply_cons   // SBind ldist_eca2011 ldist_eca2009 ldist lf_cons lcost_comply_cons // lf_cons ldist lcost_comply_cons    // has to be consistent with below!!!
	* used for original draft
	*local yvars avg_speed_noeca2011 // avg_speed_noeca2009 avg_speed_eca2009 avg_speed_eca2011  //  dist_eca2009    f_eca09_cons dist f_cons // f_eca09_cons f_cons f_eca11_cons // f_cons dist dist_eca2009 dist_eca2011 
	local rte_fe route_id
	
	if "`rte_opt'"=="all_coast" {
		local rte_cond 	`"inlist(port_agg,"SoCal","NoCal") & route_type=="Coastal" "' 
	}
				
	if "`rte_opt'"=="all_exits" {
		local rte_cond 	`"inlist(port_agg,"SoCal","NoCal") & route_type=="EnterExit" "' 
		// local rte_fe route_agg // collapsing exit routes
	}
	
	if "`rte_opt'"=="socal_coast" {
		local rte_cond 	`"port_agg=="SoCal" & route_type=="Coastal" "' 
	}
	
	if "`rte_opt'"=="socal_exits" {
		local rte_cond 	`"port_agg=="SoCal" & route_type=="EnterExit"  "' 
		local rte_fe route_agg // collapsing exit routes
	}
	
	if "`rte_opt'"=="nocal_coast" {
		local rte_cond 	`"port_agg=="NoCal" & route_type=="Coastal" "' 
	}
			
	if "`rte_opt'"=="nocal_exits" {
		local rte_cond 	`"port_agg=="NoCal" & route_type=="EnterExit" "' 
		local rte_fe route_agg // collapsing exit routes
	}

	if "`rte_opt'"=="other_coast" {
		local rte_cond 	`"port_agg=="OtherWC" & route_type=="Coastal" "' 
	}

	if "`rte_opt'"=="other_exits" {
		local rte_cond 	`"port_agg=="OtherWC" & route_type=="EnterExit" "' 
		local rte_fe route_agg // collapsing exit routes
	}
		
	capture drop tmp_samp 
	gen tmp_samp = sample_RD_regs & bad_month~=1 & vesseltype_regstr=="`vestype'" & `rte_cond'  & eca_ind<=2
	* dropping 6/2009 using bad_month
	
	
	*capture drop tmp_samp_unex
	*gen tmp_samp_unex = sample_RD_regs & vesseltype_regstr=="`vestype'" & exposed==0
	
	foreach fe_opt of local fe_opts {
		if "`fe_opt'"=="rtebyves" {
			local fe `rte_fe'#ves_id
		}
		if "`fe_opt'"=="rte" {
			local fe `rte_fe'
		}
				
		preserve
		*capture {
			*gen cut = (`run'>=0)
			*reghdfe `yvar'  `controls_resid' if tmp_samp==1 , absorb(`fe') resid( resid )		
				
			local xlabs `" 1 "2009m1" 12 "2010m1" 24 "2011m1" 36 "2012m1" 42 "2012m7" "'	
			*local xlabs_opts "format(%tqMY) angle(vertical)"
			local xlabs_opts "format(%tqMY)"
				
			* speed (in and out)
			local yscale "-10(5)5"
			reghdfe avg_speed_eca2009 ib`omit_bin'.`bin_var' `controls_resid' if tmp_samp==1 , /// 
			coefleg absorb(`fe') vce(cluster ves_id) 
			estimates store speed_in  

			reghdfe avg_speed_noeca2009 ib`omit_bin'.`bin_var' `controls_resid' if tmp_samp==1 , /// 
			coefleg absorb(`fe') vce(cluster ves_id) 
			estimates store speed_out
			
			coefplot (speed_in , `opts1' ) ( speed_out , `opts2' ) , keep(*.`bin_var')  baselevels vert yline(0, lcolor(black)) recast(scatter) ///
			xlabel(  `xlabs'  ,  `xlab_opts')  ///  
			xline( `xline_opts' ) name(speed_`rte_opt', replace) ytitle("Speed (km/hr)") ylab(`yscale' , nogrid) ///
			legend(order(2 "Within ECA" 4 "Outside ECA") ) 
			graph export "rd_plots\rd_plot_full_speed_`rte_opt'.eps", replace
			
			* distance (in eca and total)
			local yscale "-350(100)150"
			*reghdfe `yvar'  `controls_resid' if tmp_samp==1 , absorb(`fe') resid( resid )			
			reghdfe dist_eca2009 ib`omit_bin'.`bin_var' `controls_resid' if tmp_samp==1 , /// 
			coefleg absorb(`fe') vce(cluster ves_id) 
			estimates store dist_in  

			reghdfe dist ib`omit_bin'.`bin_var' `controls_resid' if tmp_samp==1 , /// 
			coefleg absorb(`fe') vce(cluster ves_id) 
			estimates store dist_out
			
			coefplot ( dist_in , `opts1' ) ( dist_out , `opts2' ), keep(*.`bin_var') baselevels vert xlabel( `xlabs' , `xlab_opts' ) yline(0, lcolor(black)) recast(scatter) ///
			xline( `xline_opts' ) name(dist_`rte_opt', replace) ytitle("Distance (km)") ylab(`yscale' , nogrid) ///
			legend(order(2 "Within ECA" 4 "Total") )
			graph export "rd_plots\rd_plot_full_dist_`rte_opt'.eps", replace

			
			* fuel
			local yscale "-50(25)25"		
			reghdfe f_eca09_cons  ib`omit_bin'.`bin_var' `controls_resid' if tmp_samp==1 , /// 
			coefleg absorb(`fe') vce(cluster ves_id) 
			estimates store fuel_in  

			reghdfe f_cons ib`omit_bin'.`bin_var' `controls_resid' if tmp_samp==1 , /// 
			coefleg absorb(`fe') vce(cluster ves_id) 
			estimates store fuel_out
			
			coefplot ( fuel_in , `opts1' ) ( fuel_out , `opts2' ), keep(*.`bin_var') baselevels vert xlabel( `xlabs' , `xlab_opts' ) yline(0, lcolor(black)) recast(scatter) ///
			xline( `xline_opts' ) name(fuel_`rte_opt', replace) ytitle("Fuel (t)") ylab(`yscale' , nogrid) ///
			legend(order(2 "Within ECA" 4 "Total") )
			
			graph export "rd_plots\rd_plot_full_fuel_`rte_opt'.eps", replace

			
		/*	
			asdsad
			
			
			*reg resid `run'
			reg resid c.`run'##i.cut  if tmp_samp==1 & ( `run'<=`bw' & `run'>=-`bw' ) 
			predict rd_line , xb
			
			reg resid c.`run'##i.cut  if tmp_samp==1 & ( `run'<=60 & `run'>=-60 ) 
			predict rd_line60 , xb
			
			reg resid c.`run'##i.cut  if tmp_samp==1 & ( `run'<=90 & `run'>=-90 ) 
			predict rd_line90 , xb
			
			reg resid c.`run'##i.cut  if tmp_samp==1 & ( `run'<=120 & `run'>=-120 ) 
			predict rd_line120 , xb
			
			reg resid c.`run'##i.cut  if tmp_samp==1 & ( `run'<=150 & `run'>=-150 ) 
			predict rd_line150 , xb
			
		*}		
		*if _rc==0 {
		
			* this is daily average
			local run_plot `run'
			collapse (first) rd_line* cut (mean) resid if tmp_samp==1 & ( `run'<=`bw' & `run'>=-`bw' ) , by(`run_plot') 
			twoway (scatter resid `run_plot' , msize(tiny)) 
			
			* this is full scatter
			*keep if tmp_samp==1 & ( `run'<=`bw' & `run'>=-`bw' )
			*local run_plot `run'
			
			egen tag = tag(`run_plot') 
			twoway (scatter resid `run_plot'  , msize(vtiny) ) ///
				(line rd_line `run_plot' if tag & `run_plot'<0 , lcolor(dkorange) ) (line rd_line `run_plot' if tag & `run_plot'>=0 , lcolor(dkorange) ) ///
					/// (line rd_line60 `run_plot' if `run_plot'<0 & `run_plot'>-60 , lcolor(blue) ) (line rd_line60 `run_plot' if `run_plot'>=0 & `run_plot'<60, lcolor(blue) ) ///
					(line rd_line90  `run_plot' if tag & `run_plot'<0 & `run_plot'>-90  , lcolor(red) )    (line rd_line90  `run_plot' if tag & `run_plot'>=0 & `run_plot'<90, lcolor(red) ) ///
					(line rd_line150 `run_plot' if tag & `run_plot'<0 & `run_plot'>-150 , lcolor(dknavy) ) (line rd_line150 `run_plot' if tag & `run_plot'>=0 & `run_plot'<150, lcolor(dknavy) ) ///
					if `run_plot'>=-`bw' & `run_plot'<`bw' , xline(0 , lpattern(dash) lcolor(black) ) ///  
					xtitle( Days from ECA )  legend(off) ytitle("`lab'")  ylab(`yscale' ,nogrid) ///
					name(`yvar'_`vestype', replace) 
			graph export "rd_plots\rd_scat_`yvar'_`rte_opt'_`fe_opt'_2009.eps", replace
		*}
		
		*/
		
		restore
		
		/*
		* around 2009 ECA
		local xline_opts_m = "xline(`=eca1my' , lpattern(dash) lcolor(black) ) " 
		reghdfe `yvar' ib592.my `controls_resid' if tmp_samp==1 & ( t_eca09<=`bw' & t_eca09>=-`bw' ) , /// 
			coefleg absorb(`fe') vce(cluster ves_id) 
		estimates store exposed
		coefplot , baselevels vert	keep(*.my)  xlabel(, angle(vertical)) yline(0, lcolor(black)) recast(scatter) ///
		xline(6 , lcolor(black) lpattern(dash)  ) name(`yvar'2009, replace) ytitle("`lab'") ylab(`yscale' ,nogrid) 
		graph export "`folder'\rd_resid_`yvar'_`rte_opt'_`fe_opt'_2009.eps" , replace
				
		* around 2011 ECA
		* only want southern california
		local xline_opts_m = "xline(`=eca1my' , lpattern(dash) lcolor(black) ) " 
		reghdfe `yvar' ib622.my `controls_resid' if tmp_samp==1 & ( t_eca11<=`bw' & t_eca11>=-`bw' ) , /// 
			coefleg absorb(`fe') vce(cluster ves_id) 
		estimates store exposed11
		coefplot , baselevels vert	keep(*.my)  xlabel(, angle(vertical)) yline(0, lcolor(black)) recast(scatter) ///
		xline(7 , lcolor(black) lpattern(dash)  ) name(`yvar'2011, replace) ytitle("`lab'") ylab(`yscale' ,nogrid) 
		graph export "`folder'\rd_resid_`yvar'_`rte_opt'_`fe_opt'_2011.eps" , replace
		*/						
									

	} 
	
}

		*foreach yvar of local yvars {
			*local lab: variable label `yvar'
			/*
			if "`yvar'"=="dist" {
				local yscale "-30(30)60"
			}

			if "`yvar'"=="f_cons" {
				local yscale "-40(20)40"
			}
			
			if "`yvar'"=="f_eca09_cons" {
				local yscale "-50(25)50"
			}
			
			if "`yvar'"=="f_eca11_cons" {
				local yscale "-40(5)10"
			}
			
			if "`yvar'"=="dist_eca2009" | "`yvar'"=="dist_eca2011" {
				local yscale "-350(100)350"
			}
			
			if "`yvar'"=="avg_speed_eca2009" | "`yvar'"=="avg_speed_eca2011" {
				local yscale "-10(5)10"
			}

			if "`yvar'"=="avg_speed_noeca2009" | "`yvar'"=="avg_speed_noeca2011" {
				local yscale "-10(5)10"
			}
			*/

/*
* around 2009 ECA (with unexposed)
reghdfe `yvar' ib592.my `controls_resid' if tmp_samp_unex==1 & ( t_eca09<=`bw' & t_eca09>=-`bw' ) , /// 
	coefleg absorb(`fe') vce(cluster ves_id) 
estimates store unexposed
coefplot exposed unexposed , vert	keep(*.my)  xlabel(, angle(vertical)) yline(0, lcolor(black)) recast(scatter) ///
xline(5 , lcolor(black) lpattern(dash)  ) name(`yvar'2009both, replace) ytitle("`lab'") ylab(`yscale' ,nogrid) 
*graph export "`folder'\rd_resid_`yvar'_`rte_opt'_`fe_opt'_2009.eps" , replace
*/

/*
* around 2011 ECA (with unexposed)
reghdfe `yvar' ib622.my `controls_resid' if tmp_samp_unex==1 & ( t_eca11<=`bw' & t_eca11>=-`bw' ) , /// 
	coefleg absorb(`fe') vce(cluster ves_id) 
estimates store unexposed11
coefplot exposed11 unexposed11 , vert	keep(*.my)  xlabel(, angle(vertical)) yline(0, lcolor(black)) recast(scatter) ///
xline(5 , lcolor(black) lpattern(dash)  ) name(`yvar'2011both, replace) ytitle("`lab'") ylab(`yscale' ,nogrid) 
*graph export "`folder'\rd_resid_`yvar'_`rte_opt'_`fe_opt'_2009.eps" , replace
*/	

/*
* dummy ECA 1 
reghdfe `yvar' ib605.my `controls_resid' if tmp_samp==1 & ( t_eca_d1<=`bw' & t_eca_d1>=-`bw' ) , /// 
	coefleg absorb(`fe') vce(cluster ves_id) 
estimates store exposed_d1
coefplot , baselevels vert	keep(*.my)  xlabel(, angle(vertical)) yline(0, lcolor(black)) recast(scatter) ///
xline(7 , lcolor(black) lpattern(dash)  ) name(`yvar'_d1, replace) ytitle("`lab'") ylab(`yscale' ,nogrid)  
graph export "`folder'\rd_resid_`yvar'_`rte_opt'_`fe_opt'_d1.eps" , replace

* dummy ECA 2 
reghdfe `yvar' ib610.my `controls_resid' if tmp_samp==1 & ( t_eca_d2<=`bw' & t_eca_d2>=-`bw' ) , /// 
	coefleg absorb(`fe') vce(cluster ves_id) 
estimates store exposed_d1
coefplot , baselevels vert	keep(*.my)  xlabel(, angle(vertical)) yline(0, lcolor(black)) recast(scatter) ///
xline(7 , lcolor(black) lpattern(dash)  ) name(`yvar'_d2, replace) ytitle("`lab'") ylab(`yscale' ,nogrid) 
graph export "`folder'\rd_resid_`yvar'_`rte_opt'_`fe_opt'_d2.eps" , replace
*/
